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Abstra ct 

A design of experiments approach has been 
implemented using computational hypervelocity impact 
simulations to determine the most effective place to add 
mass to an existing metallic Thermal Protection System 
(TPS) to improve hypervelocity impact protection. 
Simulations were performed using axisymmetric 
models in CTH, a shock-physics code developed by 
Sandia National Laboratories, and validated by 
comparison with existing test data. The axisymmetric 
models were then used in a statistical sensitivity 
analysis to determine the influence of five design 
parameters on degree of hypervelocity particle 
dispersion. Several damage metrics were identified and 
evaluated. Damage metrics related to the extent of 
substructure damage were seen to produce misleading 
results, however damage metrics related to the degree 
of dispersion of the hypervelocity particle produced 
results that corresponded to physical intuition. Based 
on analysis of variance results it was concluded that the 
most effective way to increase hypervelocity impact 
resistance is to increase the thickness of the outer foil 
layer. Increasing the spacing between the outer surface 
and the substructure is also very effective at increasing 
dispersion. 
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Introduction 

Hypervelocity impact resistance of metallic 
thermal protection systems (TPS) proposed for future 
reusable launch vehicles is an important design 
consideration. The impact velocity between orbiting 
particles (space debris and micrometeoroids) and 
orbiting vehicles is typically in excess of material sound 
speed, resulting in the term hypervelocity impact. The 
average relative velocity between space debris and 
orbiting vehicles is typically 10 km/s, but can be as 
high as 16 km/s. At these velocities, even small 
particles have sufficient momentum to significantly 
damage a spacecraft. Metallic TPS may cover a large 
percentage of the vehicle exterior surface, resulting in 
increased probability over the design life of the vehicle 
of impact with space debris and micrometeoroid 
particles large enough to penetrate the TPS. 
Penetration of the TPS may decrease the level of 
thermal protection provided on reentry, and the impact 
debris may also damage components underneath the 
TPS. Clearly, it is important to understand the impact 
resistance of metallic TPS so more impact resistant, yet 
efficient designs can be developed. 

A review of metallic TPS development evolution is 
included in Reference 1. Several of the more recently 
developed metallic TPS concepts include a foil-gage 
honeycomb sandwich on the outer surface and a thin 
foil on the inner surface that encapsulates fibrous 
insulation. One such concept, a superalloy honeycomb 
TPS, is illustrated in Figure L As shown in the figure, 
it consists of a foil-gage honeycomb -sandwich, several 
inches of fibrous insulation, an inner foil backing, and 
the substructure. The honeycomb sandwich is made of 
Inconel 617 for vehicle regions that experience high 
temperatures during reentry and titanium for low 
temperature regions. The inner foil backing is made of 
titanium for all vehicle regions. The hypervelocity 
impact analysis of TPS with Inconel 617 honeycomb- 
sandwich and inner titanium foil is the focus of this 
paper. 

The configuration of the metallic TPS design in the 
center of a panel, away from the ed^es and attachments, 
is very similar to a Whipple Shield 2 , a structure used to 
protect satellites from space debris. As shown in Figure 
2, the basic Whipple Shield design consists of one or 
more thin walls separated by a standoff distance from 
the structure being protecting. The center of a TPS 
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panel (excluding the fibrous insulation and honeycomb 
core) is similar to a Whipple Shield, consisting of three 
thin foil layers with spacing between them. Therefore, 
the TPS panel may already be in a good configuration 
to shield against orbital debris. But, with modification 
of key design parameters it may be possible to 
significantly enhance impact resistance. The objective 
of this study is to identify the most mass efficient way 
to improve impact resistance. 

Method 

In this study a computational model of a simplified 
configuration, representative of the center of a metallic 
TPS panel, was developed. Damage predicted using the 
model was compared to experimental results. The 
mode! was then used in a statistical sensitivity analysis 
to quantify and rank influences of key design 
parameters on impact resistance. 

Model Description 

Hypervelocity impact models were generated in 
CTH 3 , a shock-physics code developed by Sandia 
National Laboratories. CTH can be used to solve 
problems involving multiple materials, large 

deformations and strong shocks. A finite volume 
technique is implemented with Eulerian mesh to solve 
the conservation of mass, momentum, and energy 
equations. Only mechanical effects are considered in 
this study. Electromagnetic, chemical, and heat 

conduction effects are neglected. 

Key features of the model developed for the 
current study are summarized in Figure 3. An 

axisymmetric model was generated to reduce 

computational time and memory requirements, allowing 
workstation-class computer resources to be used. The 
honeycomb -sandwich -panel outer facesheet and inner 
facesheet will be referred to as Layer 1 and Layer 2, 
respectively. The honeycomb core between Layer 1 
and Layer 2 was omitted. Previous computational 
studies show that impact normal to the outer surface at 
different locations relative to the core can produce 
significantly different responses, but that omitting the 
core produces an average response. 4 The Saffil fibrous 
insulation layer was also neglected because it could not 
be readily modeled with CTH. The low density 
insulation would attenuate the effects of the impacting 
particle somewhat, so that neglecting the insulation will 
produce a conservative result. The titanium foil layer 
and aluminum substructure were modeled, and will be 
referred to as Layer 3 and substructure, respectively. 

Material modeling in hydrocodes can be broken 
into three categories: volumetric response, resistance to 
distortion, and failure. 5 Material volumetric response, 
or resistance to compression, is predicted using an 
equation-of-state (EOS). Of the three material 


modeling categories, volumetric response has the 
largest influence on hypervelocity impacts. 6 For this 
reason, highly detailed tabularized EOS are used with 
CTH that report pressure and energy as a function of 
density and temperature. The EOS used are valid 
through the range of densities and temperatures 
experienced in a hypervelocity impact, and account for 
multiple material phase changes. Material resistance to 
distortion, or constitutive response, is known to be in 
general a second order effect at the pressures involved 
in hypervelocity impacts. 6 This is because the impact 
pressure is many times in excess of material strength. 
For this reason, a simple elastic -plastic model with Von 
Mises yield surfaces was used for the hypervelocity 
projectile and outer TPS layers. The aluminum 
substructure was modeled using the more detailed 
Johnson-Cook viscoplastic model, since it was expected 
that impact pressure would be much Lower in the 
substructure. Failure in tension was modeled using the 
standard CTH option in which tension is relieved by 
introducing a void into computational cells when the 
hydrostatic tensile stress exceeds fracture strength. 7 

Table 1 contains the EOS, yield strength and 
fracture strength values used. The SESAME 8 
aluminum EOS was used for both A1 1100-0 and A1 
2024-T81, with a density correction factor used with AI 
2024-T8I to account for the lower density due to 
alloying. A tabularized version of the Analytical 
Equation Of State (ANEOS) 9 was used for titanium foil 
layers. The EOS for Inconel 617 was approximated by 
using the EOS for nickel with a density correction 
factor of 1 .06 to correct for the density difference. 

Materials in Table 1 with a numerical entry under 
yield strength used an elastic -plastic model with Von 
Mises yield surface for material strength. The 
substructure was modeled using the Johnson-Cook 
viscoplastic model 10 with default values for 2024-T35I 
aluminum, since no constitutive model was available 
for 2024-T81 . Failure in tension was modeled using the 
standard CTH option in which tension is relieved by 
introducing a void into computational cells when the 
mean principal stress exceeds the fracture strength. The 
fracture strength for 2024 aluminum is estimated from 
values for 2219 aluminum given by Hertel. 2 Fracture 
strength for titanium is derived from the dynamic flow 
stress reported by Babcock 11 multiplied by a factor of 
two, as suggested by Silling. 7 Values for the dynamic 
flow stress of Inconel 617 were not available. 
However, dynamic yield strength % typically on the 
order of twice the static yield strength. Fracture 
strength for Inconel 617 was obtained by multiplying 
the yield strength 12 by a factor of two (to approximate 
dynamic yield strength) and then multiplying by an 
additional factor of two (as suggested by Silling 7 ). 

Axisymmetric models were generated in which a 
rectangular mesh of computational cells was created in 
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a radial plane. Thus, each rectangular computational 
cell represents a toroid about the axis of symmetry. As 
shown in Figure 4, a centerline symmetric boundary 
condition was imposed along the axis of symmetry, 
preventing mass, momentum, stress deviator, and 
energy fluxes through the axis of symmetry. Absorbing 
boundary conditions were used along all other 
boundaries, allowing material to flow out of the 
computational domain. Material layers were un- 
restrained. Since impact occurs at velocities in excess 
of material sound speeds, the conditions at the edges of 
a material layer will not influence the response during 
the initial microseconds of impact, the time frame of 
interest. 

It is standard practice to use four or more 
computational cells through the thickness of a layer in 
order to allow adequate resolution of the shock 
structure. 13 In addition, to insure an accurate numerical 
solution computational cells should be square in 
impacted areas and cell dimensions should only be 
varied gradually. A two-step Eulerian solution scheme 
is used in CTH which requires that both material and 
void be meshed. The geometry of the model in this 
study, with large spacing between relatively thin layers, 
requires an extremely dense mesh. In order to make 
solution of the model on the computer resources 
available more practical, it was decided to only require 
two computational cells through the wall thickness, the 
rationale being that the shock structure in the thin layers 
may not be captured completely, but that the shock 
structure generated in the larger, more massive 
incoming particle will be generated with reasonable 
accuracy. In a previous study, it was shown that this 
assumption was reasonable for problems with this 
geometry. 14 This reduces the required number of 
computational cells by a factor of 4 and total runtime by 
a factor of 5. Runtime is on the order of 1 CPU-day for 
this modeling scheme, running on an SGI MIPS 
R 10000 CPU. 

During a hypervelocity impact simulation, a large 
amount of low density gas is generated, which tends to 
dramatically increase computational run time and may 
lead to numerical difficulties. To avoid computational 
problems, a material filter was used to remove material 
with both density below 0.005 gm/cm 3 and pressure 
below 0.1 MPa. Discarding material of such low 
density is assumed to have a negligible effect on the 
simulation. 

Comparison with Experimental Results 

Figure 5 shows an impact simulation of an all- 
Titanium TPS panel, with images at several times 
superimposed. Impact was with a 3/16 in. particle with 
a 7.1 km/s velocity for correlation with experimental 
data reported in Reference 15. Experimental 
parameters are listed in Table 2. The temperatures 


indicated in Figure 5 are the melting point (930 K) and 
vaporization temperature (2300 K) of aluminum at 
atmospheric pressure. This is to allow the state of the 
material to be estimated as the impact progresses. It 
can be seen that the outer foil layer, despite being much 
thinner than the projectile diameter, has shocked the 
projectile at 1 .0 fisec, causing spalling to occur. 

Spalling is a tensile failure resulting from the reflection 
of compressive waves at free surfaces. The debris 
cloud is composed mostly of solid particles at this 
point. Upon impact with the second layer, the debris 
cloud is further shocked. At this time the leading edge 
of the debris cloud is composed of vaporous debris, 
w ith solid debris located at the center of the debris 
cloud. By 8.5 jisec, the debris cloud has dispersed 
considerably and impact is initiated with Layer 3. 
Predicted substructure damage from the computational 
models was compared to experimental data. As shown 
in Figure 6, agreement was very good for the titanium 
TPS model The picture and computational plot are at 
the same scale, however the picture is a top down view 
of the substructure and the computational plot is a side 
view. The experimental hole diameter w 7 as 0.94 cm, 
while hole size obtained computationally was 0.78 cm. 
The Inconel 617 model significantly over-predicted the 
substructure hole size, as shown in Figure 7. The most 
likely explanation for the discrepancy is that the 
adjusted nickel EOS, used to represent Inconel 617, did 
not adequately represent the behavior of Inconel 617. 

In addition, several assumptions were required to obtain 
a fracture strength for Inconel 617, which may have 
resulted in inaccuracies. It is also possible that, in the 
experimental data, the particular impact location 
relative to the foil walls in the honeycomb core 
increased dispersion of the impacting particle, resulting 
in less substructure damage. The good correlation with 
experiment for the titanium honeycomb sandwich gives 
confidence that the modeling approach is valid. 

Although the lack of accurate material properties for 
Inconel 617 may effect the predicted damage level, it is 
assumed that the Inconel model would produce results 
that show the correct trends. 

Statistical Sensitivity Analysis. 

The primary objective of this study is to determine 
the most effective place to add mass to an existing 
metallic TPS to improve hypervelocity impact 
protection. The effects of a large number of parameters 
or factors needed to be quantified and compared. An 
efficient method was needed to assess the relative 
importance of many factors. The conventional 
approach would be to investigate one factor at a time, 
with the remaining factors held constant. This method 
is acceptable for a few variables that do not change 
much from their baseline values. Moreover, this “one- 
factor-at-a-time” method assumes that the effect of each 
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variable is independent of variations of the other 
variables. 

A factorial design approach can be used in which a 
few discrete values are chosen for each design 
parameter and all possible combinations are evaluated. 
However, the number of analyses required can quickly 
become unmanageable. 

The method selected for this study borrows 
techniques commonly used in robust design. 1617 
Design parameters (independent variables) a-e varied 
simultaneously in a matrix of experiments (or in this 
case, analyses) defined by orthogonal arrays. The 
effects of a large number of design parameters, 
evaluated over the range of interest, can be determined 
efficiently using this method. The method used in this 
study follows closely the procedures presented in a 
chapter entitled “Matrix Experiments Using Orthogonal 
Arrays” in Reference 17 as well as similar procedures 
in Reference 16. 

The first step in statistical sensitivity analysis is to 
identify the design parameters of interest. Care must be 
taken to select a meaningful and reasonable range of 
variation for each design parameter. The range of 
variation can affect how much each factor influences 
the results. The selection of the factors and their ranges 
usually depends on engineering experience and 
judgment, but a good selection is essential to formulate 
the problem properly for meaningful results. Several 
levels, usually 2 to 4, are selected to cover the range of 
each design parameter. 

The next step in statistical sensitivity analysis is to 
determine what experiments or analyses need to be 
performed to determine the effects of the selected 
factors on the performance of the system. Orthogonal 
arrays are used to define the number of analyses to be 
performed and the combination of design parameter 
levels for each analysis. A limited number of standard 
orthogonal arrays 16 are available to accommodate 
specific numbers of design parameters with various 
levels per design parameters. In the orthogonal array 
there is a column for each design parameter, and each 
row is a particular combination of levels for each 
design parameter. The columns of the array are 
mutually orthogonal; that is, for any pair of columns, all 
combinations of design parameter levels occur, and 
they occur an equal number of times. 17 The smallest 
standard orthogonal array is chosen to accommodate 
the number of design parameters and design parameter 
levels. Any extra columns in the orthogonal array can 
either be ignored or in some cases used to estimate the 
effects of interactions between design parameters. 

An experiment, or in this case an analysis, is 
performed for every row in the orthogonal matrix with 
the combination of design parameter levels defined in 
that row. Interpretation of the results is straightforward. 
The result or results for each analysis are tabulated and 


the overall means for all analyses are calculated. The 
means for each level of every design parameter are also 
calculated. Each column of the orthogonal matrix 
contains an equal number of experiments at each level 
of the design parameter associated with that column. 
The results associated with each level of that design 
parameter are averaged to obtain the associated means. 
The effect of a design parameter level is defined as the 
deviation it causes from the overall mean. Therefore, 
the effect of each design parameter level can be 
obtained by subtracting the overall mean from the mean 
associated with the particular level of that design 
parameter. This process of estimating the design 
parameter effects is sometimes called analysis of means 
(ANOM). 

Using ANOM the effect of each design parameter, 
called the main effects, can be determined 
independently. A simple additive or superposition 
model is assumed for the response 7]: 

T] = overall mean + ^(design parameter effect ) + e ( 1 ) 

design 
jwra meter 

where £ is the error of the additive approximation. 
Thus the response for any combination of design 
parameter levels can be estimated using this equation. 
In a typical design of experiments analysis this error 
term would also include a contribution associated with 
the error in repeatability of measuring the response, Tj, 
for a given experiment. However, for the current 
application the response is obtained from analysis, so 
error associated with measurement repeatability does 
not exist. 

The results of ANOM can provide considerable 
insight into the effect of the various design parameters; 
however, an analysis of variance (ANOVA) can give a 
more accurate indication of the relative importance of 
the design parameters and provides a means of ranking 
the factors in order of importance. ANOVA can be 
used to determine the contribution of each design 
parameter to total variation from the overall mean 
value. The sum of squares of the differences from the 
mean for all the levels of a particular design parameter 
(Column 5 of Table 4) provides a measure of how much 
that design parameter affects the result over the 
specified range. The percentage that this sum of 
squares value contributes to the total for all design 
parameters (% Total SS in Table 4) gives a measure of 
the relative importance of that particular design 
parameter. 

As part of the ANOVA, the error associated with 
the additive assumption can be estimated. The method 
used for the error estimate depends on the number of 
design parameters and design parameter levels, and the 
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orthogonal matrix used in the statistical sensitivity 
analysis. 

The sum of squares due to error 17 can be calculated 
using the following relationship: 

Sum of squares due to error = 

(grand total sum of squares) 

-(sum of squares due to mean) ( 2) 

-(sum of squares due to design parameters) 

where the grand total sum of squares is the sum of the 
squares of all values for a particular result, the sum of 
squares due to the mean is the overall mean squared 
times the number of analyses (or number of rows in the 
orthogonal matrix), and the sum of squares due to 
design parameters is the sum of squares of all the 
design parameter effects. There are degrees of freedom 
associated with each of the quantities in Equation 2. 
The degrees of freedom for the grand total sum of 
squares are the number of rows in the orthogonal 
matrix. There is one degree of freedom associated with 
the mean. Each design parameter has one less degree 
of freedom than the number of levels for that design 
parameter (d.o.f in Table 4). Therefore the degrees of 
freedom for the error can be calculated as follows; 

(degrees of freedom for error) = 

(number of rows in orthogonal matrix - 1) (3) 

- (sum of degrees of freedom for all 
design parameters) 

The degrees of freedom for the error must be 
greater than zero for Equation 2 to be useful. 

If there are zero degrees of freedom available for 
calculating error, as is the case in this study, a different 
method must be used to estimate the additive error. An 
approximate estimate of the sum of the squares due to 
error can be obtained by pooling the sum of squares 
corresponding to the design parameters having the 
lowest mean square. A rule of thumb 17 is to use the 
sum of squares corresponding to the bottom half of the 
design parameters (as defined by lower mean square) 
and the degrees of freedom corresponding to those 
design parameters. Pooled error estimates are used for 
the results shown in this paper. 

Once the sum of square due to error and the 
degrees of freedom for error have been calculated, the 
error variance can be estimated as follows: 

Error variance = (sum of squares due to error)/ 

(degrees of freedom for error) (4) 

The variance ratio, F, is a measure of how 
important the effects of a design parameter are 
compared to the error. 


F = (mean square due to a factor)/(error variance) (5) 

The mean square due to a design parameter (Mean 
SS in Table 4) is the sum of squares of the differences 
from the mean for a design parameter (SS) divided by 
the degrees of freedom associated with that design 
parameter (dof). 

Design Parameters 

In all sensitivity analyses, a particle diameter of 
3/16 in. was assumed with a normal impact velocity of 
7.1 km/s. From a series of experimental impacts on 
metallic TPS, 15 it was found that 1/8 in. particles could 
be stopped by most TPS configurations without 
modification. All tests with larger Va in. diameter 
particles resulted in penetration of the TPS and 
substructure, and it was clear from the level of damage 
that significant modification of the TPS would be 
required to provide protection against a particle of this 
size. However, substructure damage in tests with 3/16 
in. diameter particles were very sensitive to changes in 
TPS parameters, making this particle size interesting for 
study in a sensitivity analysis. 

The thicknesses of each of the three layers, as well 
as the spacing between Layer I and Layer 2 and the 
spacing between Layer 1 and Layer 3, were selected as 
design parameters. Figure 8 lists the five design 
parameters, labeled A through E along with the design 
parameter levels studied. Levels were chosen that 
would encompass a reasonable design range for 
metallic TPS panels. In addition, design parameters A 
through C were calculated so that the difference in mass 
between Level 1 and Level 2 would be the same for 
each, 0.1 1 lb per panel. This is to determine the most 
effective location in the TPS to add mass. Level 1 and 
Level 2 values for design parameters D and E were 
chosen to encompass the design range expected for 
metallic TPS panels, and do not represent equal mass 
differences in a simple way. The change in mass for 
parameter D is dependent on the density of the 
honeycomb core. For the honeycomb used in 
experimental tests 18 , a 0.28 in. increase in honeycomb 
thickness results in a 0.19 lb/ft 2 weight increase. The 
change in mass for parameter E is dependent on the 
specific metallic TPS design being used, and the 
method of varying the standoff between Layer 1 and 
Layer 3. The following discussion will focus on the 
metallic TPS panel design shown in Figure 1. If the 
insulation thickness is increased and the beaded side 
wall elongated to accommodate the 1 in. increase in 
overall TPS thickness between Level 1 and Level 2 of 
design parameter E, the weight increase, calculated 
using data from Reference 18, is 0.314 lb/ft 2 . However, 
the increase in insulation weight accounts for 83% of 
this weight, and is not required. With modification, the 
TPS panel could be constructed with a spacer 
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underneath the insulation layer that would allow TPS 
panel thickness to be varied without increasing 
insulation thickness, potentially reducing the weight 
increase related to parameter E. There are, however, 
other less obvious weight penalties related to increasing 
the TPS panel thickness. Increasing design parameter E 
from Level 1 to Level 2 significantly effects the outer 
mold line of the vehicle, so that the impact on vehicle 
weight or performance may be greater than just the TPS 
weight increase. 

Orthogonal arrays are used to define the number of 
analyses to be performed and the combination of design 
parameter levels for each analysis. A limited number of 
standard orthogonal arrays are available to 

accommodate specific numbers of design parameters 
with various levels per parameter. 16 An Li 6 two level 
orthogonal array, shown in Table 3, was selected for 
statistical sensitivity analysis since in this case it 
provides high resolution for the determination of design 
parameter sensitivity. The L| 6 array contains 15 

columns and requires 16 experiments or analyses to be 
performed. Design parameters can be assigned to 
certain columns dependent on the orthogonal array 
chosen. In this study, columns 1, 2, 4, 8, and 15 were 
used for design parameters A through E. A signal 
strength, which can be thought of as a signal to noise 
ratio, was calculated for each column to reveal the 
system sensitivity to variation of the design parameters. 
Signals were also obtained for interactions between the 
design parameters. The other columns in Table 3 
contain interactions between the design parameters, 
where “x” represents interaction. For example, AxB, 
located in column 3, represents the interaction between 
design parameters A and B. The columns containing 
interactions are confounded, meaning that multiple 
parameters are contained in the same column. In this 
case, each double interaction parameter is confounded 
with a triple interaction parameter. In the previous 
example, column 3 contained the interaction between 
design parameters C, D, and E in addition to the 
interaction between design parameter A and B. This 
indicates that predicted interaction sensitivities will be 
less accurate than predicted design parameter 
interactions. The emphasis in this research was on the 
effect of design parameter variation, and the design 
parameters were found to be more influential than the 
interactions. Each row represents a particular analysis, 
where the 44 1” and “2” entries in the remaining columns 
represent the design parameter levels. As shown in 
Table 3, there are 16 rows, trials 1 through 16, where 
each row represents a single analysis. The numbers in 
the row, either 1 or 2, represent the level of the design 
parameters for the analysis. Numbers listed in columns 
that only contain interactions are only used in 
interpretation of the results and have no influence on 
the design parameter levels analyzed. 


Damage Metrics 

Use of statistical analysis requires identification of 
a response, T| , to be measured as a “factor of goodness”, 
a way of quantifying the positive effect of varying 
design parameters. Careful selection of this factor is 
essential to obtaining meaningful results, and the 
parameter chosen should directly represent the physics 
involved in the process. 16 

In this study, TJ relates to the effectiveness of the 
TPS at protecting the substructure. It is a measure of 
the damage potential to the substructure, or in other 
words, a damage metric. Two forms of damage metrics 
were investigated: predicted substructure damage and 
degree of dispersion of the incoming particle. 

Predicted substructure damage using hole sizes was 
the more intuitive damage metric, however statistical 
analysis using substructure hale size as a factor of 
goodness produced misleading results. Figure 9 shows 
three conceptual impact cases that reveal why this is the 
case. The first case, labeled “No Bumper”, shows 
impact on an unprotected substructure, resulting in a 
hole approximately the size of the incoming particle. 
The second case, labeled “Thin Bumper”, uses a 
Whipple Shield that results in particle break-up, but 
does not sufficiently disperse the particle debris cloud. 
The result is a substructure hole size that is larger than 
if the substructure was unprotected! Finally, the third 
case, labeled “Thick Bumper” represents an effective 
Whipple Shield that disperses the debris cloud to the 
extent that the substructure withstands the impact. If 
these shield designs were ranked based on substructure 
hole size the “No Bumper” design would be ranked 
above the “Thin Bumper” design, even though the thin 
bumper design is reducing and dispersing the impact 
energy. ^ 

The other form of damage metric identified focuses 
on the degree of dispersion of the incoming particle. In 
this damage metric the state of the particle debris cloud 
after impact with TPS, prior to substructure impact, is 
used. Identification of a single number that would 
adequately describe the degree of debris cloud 
dispersion was difficult. It was decided to measure 
what will be referred to as “Percent Momentum 
Radius” for several different percentage values. 
Percent momentum radius is defined as the radius, 
measured from the axis of symmetry, that contains a 
specified percentage of the initial projectile axial 
momentum. This concept is illustrated in Figure 10 for 
an idealized impact scenario, where 50% and 100% 
momentum radii are shown as a function of time. 
Percent momentum radius is a time dependent quantity, 
which requires that a consistent time be used for 
measurement in all of the cases examined. Because the 
geometry of the cases varied considerably, a fixed time 
was not used (for example 10 jisec.). Instead, percent 
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momentum radius was recorded at the initiation of 
substructure impact, which varied from 10 to 16 ^isec 
after initial particle impact. 

Percent momentum radius was calculated over the 
active computational cells only. To improve run-time 
and eliminate potential numerical problems, the 
computational cells containing Layer 1 and Layer 2 
were removed after impact. In addition, momentum 
was lost out of the sides of the computational domain, 
which was sufficiently dispersed to be of no threat to 
the substructure but does result in a significant drop in 
total momentum. Due to these momentum losses, the 
largest percent momentum that could be used for all 
cases was 30%. In addition, 10% and 20% momentum 
radius were calculated. 

Percent momentum radius was used in preference 
to substructure hole size beause it was more directly 
related to the physics of particle break-up. Use of this 
parameter resulted in the ranking of TPS designs based 
on their ability to disperse incoming particles, which 
had the added advantage of uncoupling the 
computational results from the specific substructure 
material and geometry simulated. 

Results 

TPS Design Guidelines 

After completing the analyses specified in Table 3, 
statistical sensitivity analysis was performed for several 
different damage metrics. Table 4 shows representative 
results for the 30% momentum radius damage metric. 
Design parameters and interactions (these are the 
columns from Table 3) are listed in order of influence 
to the damage metric. Design parameter E has the most 
effect on the 30% momentum radius damage metric. 
Column 2 in Table 4 lists a description of the design 
parameters. For example, design parameter E is the 
spacing between Layer 1 and Layer 3. Columns 3 and 
4 list the difference from mean of the damage metric at 
design parameter Level 1 and Level 2 respectively. In 
this case, design parameter E results in an average - 
0.17 cm decrease in the 30% momentum radius when 
set to 2.28 in. and a 0.17 cm increase when set to 3.28 
in.. Columns 5 and 6 report the sum of squares of the 
damage parameter and percent total of the sum of 
squares. Column 7 lists the degrees of freedom of each 
design parameter and interaction, which is 1 in all cases 
because the array is limited to two levels. Column 8 is 
the mean sum of squares, which was calculated by 
dividing the sum of squares by the number of degrees 
of freedom. 

The final step in statistical sensitivity analysis is to 
calculate the variance ratio, F. The variance ratio, 
similar to a signal to noise ratio, measures the 
importance of the effects of a design parameter 
compared to the error or noise in the results. The 


variance ratio, F, calculated using equation 5, is shown 
in column 9 of Table 4. 

Table 5 lists the two most influential design 
parameters for each of the damage metrics examined. 
Statistical sensitivity analyses using 20% and 30% 
momentum radius both predict that a l .0 in. variation of 
design parameter E will have the most influence on 
debris cloud dispersion. The second most influential 
design change was the 0.0025 in. variation of design 
parameter A. Both cases predict that increasing the 
value of these parameters will increase debris cloud 
dispersion and thereby provide better protection from 
hypervelocity impacts. 

Even though it is predicted that design parameter E 
is more influential than design parameter A, it is likely 
that from a design standpoint design parameter A is the 
better parameter to change. Change of design 
parameter E, although more influential than changing 
A, results in a 0.314 lb increase in TPS weight, while 
the increase of parameter A from Level 1 to Level 2 
results in a 0.11 lb increase in TPS weight. Although 
modifications could be made to the TPS design to 
reduce the weight increase associated with increasing 
parameter E, increasing the thickness of the TPS panel 
influences the outer mold line of the vehicle which will 
effect vehicle weight and performance. Increasing 
design parameter A, the thickness of Layer l , provides 
additional benefits to the vehicle from an operability 
standpoint. Increasing the outer facesheet thickness 
increases the ability of the vehicle to fly through rain, 
resist ground hail strikes, and resist damage during 
assembly. 

Further examination of Table 5 reveals that there is 
a discrepancy between predictions using 10% 
momentum radius as a damage metric and those using 
larger radii, such as 20% momentum radius and 30% 
momentum radius. Although there is agreement that 
design parameter E is influential, using the 10% 
momentum radius damage metric results in the 
prediction that a 0.28 in. variation in design parameter 
D is the most influential, and that the value of the 
parameter should be decreased to improve impact 
protection. This non-intuitive result led to a hypothesis 
that there must be some phenomena leading to the 
difference of predictions resulting from using 10% 
momentum radius as a damage metric, as opposed to 
the 20 or 30% momentum radius, and that the 
difference was probably related to the fact that 10% 
momentum radius measurements focus on a narrower 
region of the impact debris cloud than 20% and 30% 
momentum radius measurements. Analysis results 
were examined to obtain a better understanding of the 
phenomena leading to this discrepancy, and are 
discussed in the following section. 
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Modeling Artifacts 

Examination of density and momentum data during 
the penetration of Layer 2 revealed a likely source for 
the discrepancy between the 10% damage metric, which 
focuses on a narrow region of the computational 
domain, and the 20% and 30 % damage metrics, which 
focus on a larger area. Momentum readings along the 
axis of symmetry were found to gradually decrease 
after impact with Layer 1 , as expected, but after impact 
with Layer 2 momentum along the centerline increased 
rapidly. This is shown in Figure 11, which plots 

annular momentum at 0.0, 3.0, and 5.0 jisec. Annular 
momentum is defined in Figure 12, and is: 

n 

Annular Momentunij) = lW,j)IA( 0 (6) 

M 

where i and j are the radial and axial computational cell 
indices, L represents the axial momentum of the 
specified cell, and A represents the cross-sectional area 
of the annulus formed by the computational cell at the 
specified radial index. At T=0.0 jLisec impact is 
initiated with Layer I. At T=3.0 psec penetration of 
Layer 2 is in progress, and the annular momentum has 
decreased substantially. However, at T=5.0 jLisec, after 
penetration of Layer 2, a spike is seen in the annular 
momentum at radius of 0.0 cm, along the axis of 
symmetry. Examination of debris cloud radial velocity 
during penetration of Layer 2 revealed that as Layer 2 
deforms it imparts a significant negative radial velocity 
to portions of the debris cloud, as shown in Figure 13. 
This “funnel effect” directs portions of the debris cloud 
to the axis of symmetry. In the axisymmetric model 
being used, the radial velocity of the funneled debris is 
set to zero when it reaches the axis of symmetry, 
causing debris to accumulate along the axis. The funnel 
effect combined with the centerline mass accumulation 
results in the increase of momentum along the 
centerline. 

Estimation of the upper bound value of centerline 
mass accumulation revealed that from 0 to 4% of the 
initial axial momentum was concentrated along the 
centerline, depending on the configuration of Layer 1 , 
Layer 2, and the spacing between these layers. In order 
to determine if this was areal phenomena or an artifact 
of the modeling technique a three-dimensional (3D) 
analysis of Layer 1 and Layer 2 was created. Results 
from the 3D model predicted localized funnel effects 
but did not predict the overall funnel effect about the 
centerline of the impact shown in Figure 13. Mass 
accumulations occurred in the 3D analysis, but were not 
limited to the centerline and quickly dissipated. Figure 
14 shows the axial momentum distribution at T=0.0, 
2.0, 4.0, and 7.0 jusec., which corresponds to the debris 
cloud state at the initiation of impact, after impact with 
Layer 1, during impact with Layer 2, and 3 jisec after 


impact with Layer 2, respectively. The axial 
momentum distribution plotted is the summation of 
computational cell axial momentum through ,y 
(corresponding to axial direction in axisymmetric 
models) at a fixed x, z location. The total momentum is 
then divided by the cell area in the x-z plane. As shown 
in Figure 14, debris cloud momentum is dispersed after 
impact with Layer 1. Immediately after impact with 
Layer 2, at T= 4.0 jisec, the axial momentum increased 
at several locations, which is probably due to the funnel 
effect. However, as shown in the axial momentum at 
T=7.0 jisec, the momentum increase is temporary, and 
the mass accumulation quickly dissipates. So, the 3D 
model verified that funnel effects lead to mass 
accumulations, however mass accumulations are not 
limited to the impact centerline, are not as large as in 
the axisymmetric model, and quickly dissipate. 

The funnel effect and centerline mass accumulation 
explains the discrepancy between the statistical 
sensitivity analysis results using 10% momentum radius 
and larger momentum radii. Upper bound estimates of 
the axial momentum due to centerline mass 
accumulation vary from 0 to 4% of the initial axial 
momentum depending on the TPS configuration. 
Therefore, the 10% momentum radius is significantly 
influenced by the centerline mass accumulation, where 
the degree of influence will be much lower for 20% and 
30% momentum radii. 

Conclusion 

A design of experiments approach was used with 
computational hypervelocity impact analyses to 
determine the most effective place to add mass to an 
existing metallic Thermal Protection System to improve 
hypervelocity impact protection. The choice of damage 
metric, used as a “factor of goodness” to rank 
compeating designs, was critical. 

Use of percent momentum radius as a damage 
metric provided a good indication of the level of 
dispersion of the hypervelocity particle prior to impact 
with the substructure. Use of 20% or 30% momentum 
radii in sensitivity studies produced results that 
correspond to physical intuition. Based on these 
results, it is clear that the best place to add mass to 
increase hypervelocity impact protection is the 
outermost surface of the TPS panel, referred to as Layer 
1 , Adding mass to the outer surface not only improves 
hypervelocity impact protection, but increases the 
resistance of the panel to rain, hail, and handling 
damage. The spacing between the outermost layer and 
the substructure was also influential, and increasing this 
spacing will increase impact protection. However, it is 
likely that the penalty will be larger, both due to TPS 
weight directly related to the thickness of the panel, and 
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vehicle weight or performance that is also related to the 
panel thickness. 

It was found that use of substructure hole size as a 
damage metric is potentially misleading because it is 
specific to an assumed substructure design and does not 
directly indicate the level of dispersion in the debris 
cloud. This may result in marginally effective 
protection systems being ranked lower than no 
protection system at all, since the marginal systems will 
break up the hypervelocity particle, potentially 
producing larger substructure holes. 

The centerline mass accumulation observed in 
axisymmetric analysis results, termed funnel effect in 
this paper, is not realistic for this problem. Use of 3D 
models showed that the funnel effect predicted by 
axisymmetric models is possible, however the resulting 
mass accumulations are not limited to the impact 
centerline, were not as large, and were non-persistent. 
The discrepancy between 10% momentum radius 
sensitivity analysis results and the 20% and 30% 
momentum radius sensitivity results is attributed to the 
centerline mass accumulation. Use of 10% momentum 
radius as a damage metric produced unintuitive results 
since the focus of this measurement is the central 
portion of the model influenced by the centerline mass 
accumulation. It is concluded that 20% and 30% 
momentum radii measurements give a better indication 
of the degree of dispersion in the debris cloud. 

Based on these findings, it is clear that the best 
place to add mass to increase hypervelocity impact 
protection is the outermost surface of the TPS panel. 
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Table 1: CTH Material Parameters 


Material 

Equation of State 

Yield Strength 
(MPa) 

Fracture 
Strength (MPa) 

11 00-0 A1 

SESAME Aluminum 

34 

-100 

2024-T81 Al 

SESAME Aluminum 

Johnson-Cook 

-1500 

Ti 6-4 

-AN EDS Ti fi- 4 



■HHlHeim 

ANEOSNI 

35Q 

-1400 


Table 2: Experimental Parameters 



Inconel TPS 

Titanium TPS 


INC 617 

Ti-6AL-4V 

IGauae fin) 

0.005 

0.005 


INC 617 

Ti-6AL-4V 

Cell Size (in) 

0.25 


Gauge (in) 

0.002 


DeDth (in) 

0.28 

0,28 

Inner Facesheet 

IHNK3BH 


Material 

■ 

Ti-6AL-4V 

Gauae fin) 

HExIiSlHi 

0.005 

Insulation 



Material 

Saffil 

Saffil 

Thickness (in) 

2.0 

2.0 

Density flb/ft 3, i 

1.5 

1.5 

rim i ii— i 



Material 

Ti-6AL-4V 

Ti-6AL-4V 

Gauae fin) 

0.003 

0.003 

Substructure 



Material 

2024-T81 A! 

2024-T81 Al 

Thickness fin) 

J2J- 

01 
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Table 3: LI 6 Two-Level Orthogonal Array. Design parameters located on columns 1, 2, 4, 8, 15. Two way 
and three way interactions located on remaining columns. 


Column Number 

6 7 8 9 10 11 12 13 14 15 


Trial # 

A 

4. 

B 

AxB, 

CxDxE 

C 

AxC, 

BxDxE 

BxC, 

AxDxE 

DxE, 

AxBxC 

D 

AxD, 

BxCxE 

BxD, 

AxCxE 

CxE, 

AxBxD 

CxD, 

AxBxE 

BxE, 

AxCxD 

AxE, 

BxCxD 

E 

1 

1 

“T 

1 

~r 

1 

1 

1 

_ 1 

1 

1 

1 

1 

1 

1 

1 

2 

1 

i 

1 

i 

1 

1 

1 

2 

2 

2 

2 

2 

2 i 

2 

2 

3 

1 

i 

1 

2 

2 

2 

2 

1 

1 

1 

1 

2 

2 

2 

2 

4 

1 

i 

1 

2 

2 

2 

2 

2 

2 

2 

2 

1 

1 

1 

1 

5 

1 

2 

2 

1 

1 

2 

2 

1 

1 

2 

2 

1 

1 

2 

2 

6 

1 

2 

2 

1 

1 

2 

2 

2 

2 

1 

1 

2 

2 

1 

1 

7 

1 

jH 

2 

2 

2 

1 

1 

1 

1 

2 

2 

2 

2 

1 

1 

8 

1 

2 

2 

2 

2 

1 

1 

2 

2 

1 

1 

1 

1 

2 

2 

9 

2 

1 

2 

1 

2 

1 

2 

1 

2 

1 

2 

! 1 

2 

1 

2 

10 

2 

1 

2 

1 

2 

1 

2 

2 

1 

2 

1 

2 

1 

2 

1 

11 

2 

1 

2 

2 

1 

2 

1 

1 

2 

1 

2 

2 

1 

2 

1 

12 

2 

1 

2 

2 

1 

2 

1 

2 

1 

2 

1 

1 

2 

1 

2 

13 

2 

2 

1 

1 

2 

2 

1 

1 

2 

2 

1 

1 

2 

2 

1 

14 

2 

2 

1 

1 

2 

2 

1 

2 

1 

1 

2 

2 

1 

1 

2 

15 

2 

2 

1 

2 

1 

1 

2 

1 

2 

2 

1 

2 

1 

1 

2 

16 

2 

2 1 

2 

1 

1 

2 

2 

1 

1 

2 

1 

2 

2 

r 


Table 4: Analysis of Variance for 30% Momentum Radius Damage Metric 


! 30% Momentum Radius 

Difference from Mean, cm 

KfnTTiyH 

BPiWmWI 


Hull 

r — 

1 Factor 

Level 1 Level 2 

E 

A 

DxE 
A x D 
D 

CxD 

B 

*AxE 

*C 

*B x C 
‘AxB 
*C x E 
‘BxD 
*B x E 
*A x C 

Spacing Layer 1 / Layer 3 
Thickness Laver 1 
Interaction 
Interaction 

Spacing Layer 1 / Layer 2 
Interaction 
Thickness Laver 2 
Interaction 
Thickness Layer 3 
Interaction 
Interaction 
Interaction 
Interaction 
Interaction 
Interaction 

-1.74E-01 1.74E-01 

-1.13E-01 1.13E-01 

8.19E-02 -8.19E-02 

6.31E-02 -6.31E-02 

-4. 94 E -02 4.94E-02 

4.94E-02 -4.94E-02 

-3.69E-02 3.69E-02 

3.56E-02 -3.56E-02 

-3.19E-02 3.19E-02 

-3.06E-02 3.06E-02 

-2.94E-02 2.94E-02 

2.69E-02 -2.69E-02 

1.69E-02 -1.69E-02 

1.19E-02 -1.19E-02 

-6.87E-03 6.87E-03 

4.87E-01 
2.05E-01 
1.07E-01 
6.38E-02 
3.90E-02 
3.90E-02 
2.18E-02 
2.03E-02 
1.63E-02 
1.50E-02 
1.38E-02 
1 .16E-02 
4.56E-03 
2.26E-03 
7.56E-04 

46.49 

19.56 

10.25 

6.09 
3.73 
3.73 
2.08 
1.94 
1.55 
1.43 
1.32 

1.10 
0.44 
0.22 
0.07 

H 

4.87E-01 

2.05E-01 

1.07E-01 

6.38E-02 

3.90E-02 

3.90E-02 

2.18E-02 

2.03E-02 

1.63E-02 

1.50E-02 

1.38E-02 

1.16E-02 

4.56E-03 

2.26E-03 

7.56E-04 

46.06 
19.39 
10.15 
6.04 
3.69 
3.69 

2.06 

PrsnlArl F 

: rmr (*\ 


Imwsmm 


1 (8) 

imib i ii \m 


= 

Total 

1 1 .05E+00 

[ 100.00 

1 
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Table 5: Recommended change in design 
parameters determined using different damage metrics 


Rank 

10% Momentum 
Radius 

20% Momentum 
Radius 

30% Momentum 
Radius 

1 




2 


IB !1!H 

Increase A 1 
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Fastener access 
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Fastener access 
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Saffil insulation 
Inner foil layer 


Mechanical 
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Figure 1 : Metallic TPS Panel 



Space debris 
approaches spacecraft 


Debris impacts a thin 
bumper layer - 
resulting shock 
disintegrates particle 


Debris cloud 
disperses. Impact is 
spread over larger 
substructure area. 


Figure 2: Basic Whipple Shield design principle 
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M etallic TPS 
Configuration 


Axisymmetric 

Computational 

M odel 


Outer Facesheet 
h/c Core 
Inner Facesheet 



Layer 1 

Core Neglected 
Layer 2 


Saffil 

Insulation 


Titanium Foil 



Insulation 

Neglected 


Layer 3 


9.5 


mm 


Aluminum 2024-T81 
Substructure 



Substructure 


<L 


Figure 3: Comparison of Experimental Configuration (through the thickness) to Axisymmetric 

Computational Model 


14 

American Institute of Aeronautics and Astronautics 



Figure 4: Schematic representation of axi symmetric 
finite volume model created for use with CTH. 
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Figure 5: Time elapsed view of computational 
simulation of hypervelocity impact on titanium TPS 
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Figure 6: Comparison of Predicted and Measured 
Substructure Damage for Titanium TPS 
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Figure 8: Definition of parameter levels used in sensitivity analysis 
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Figure 9: Schematic representation of hypervelocity impact with three different shielding configurations. 



Figure 10: Explanation of percent momentum radius using 
idealized diagram of axisymmetric impact problem, with 
debris cloud shown at several arbitrary times. 
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1 .OE+05 



Radius (cm) 

Figure 11: Annular momentum distribution in typical 
axisymmetric simulation, revealing centerline mass 
accumulation after impact with Layer 2 (at 5 psec.) 



Annular Momentum(i) = ^ £(/, j) I A(i) 

i 


Figure 12: Schematic showing the summation of axial 
momentum in annular elements, where L(i, j) 
represents the axial momentum of annular element 
at radial indice i and axial indice j. 
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Figure 13: Initiation of funnel effect during penetration of Layer 2 



number 




Figure 14: Axial momentum distribution at four timesteps for three-dimensional model. Momentum spikes 
at 4. psec. reveal possible funnel effects. However, funnel effects are non-persistent and do not lead to 
centerline mass accumulation as predicted in axlsymmetric model. 
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